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Abstract 

Rather  than  investigate  a particular  combustor,  a whole  class 
of  combustion  systems,  susceptible  to  damage  from  combustion 
instabilities,  is  considered-  Under  some  simple  and  realistic  as- 
sumptions (pressure  waves  reflected  from  the  combustor  bound- 
aries smaller  than  incoming  waves,  flame  stable  in  itself,  limited 
bandwith  flame  response),  it  is  demonstrated  that  a finite  dimen- 
sional approximation  to  the  open-loop  transfer  function  of  such  a 
combustion  system  satisfies  some  general  properties  (stable  zeros, 
small  relative  degree)  that  are  exploited  to  design  adaptive  ac- 
tive controllers  guaranteed  to  stabilise  the  self-excited  combustion 
oscillations.  In  particular,  for  the  practical  case  of  a combustion 
system  with  time  delay,  a completely  new  and  simple  adaptive  con- 
trol design  is  presented  and  a formal  proof  for  stability  is  given. 
The  performance  of  such  stable  adaptive  controllers  is  illustrated 
in  a simulation. 

1.  Introduction 

In  order  to  meet  stringent  emission  requirements,  combustors 
are  increasingly  being  designed  to  operate  in  a lean  premixed 
mode.  Although  this  reduces  the  NOx  emissions,  it  has  the  disad- 
vantage that  premixed  flames  are  particularly  susceptible  to  self- 
excited  oscillations,  and  the  associated  large-scale  prassurc  waves 
can  cause  structural  damage.  Active  control  provides  a way  of 
extending  the  stable  operating  range  of  a combustion  system,  by 
interrupting  the  damaging  interaction  between  acoustic  waves  and 
unsteady  combustion.  However,  to  be  useful  in  practice,  an  active 
controller  needs  to  be  effective  across  a range  of  operating  con- 
ditions. An  efficient  approach  is  to  use  an  adaptive  controller  in 
which  the  controller  transfer  function  is  continually  altered  as  the 
engine  condition  changes. 

They  are  already  some  algorithms  that  describe  how  to  update 
the  controller  parameters.  The  most  popular  adaptive  schemes 
used  for  active  control  of  combustion  instability  is  the  Least  Mean 
Squares  (LMS)  algorithm  applied  to  an  IIR  (Infinite  Impulse  Re- 
sponse) filter  [5], [17], [12], [11].  The  LMS  is  very  attractive  because 
it  does  not  require  any  theoretical  model:  the  combustion  process 
is  considered  as  a ‘black  box’  and  is  learnt  during  a system  iden- 
tification procedure,  performed  off-line  [5]  or  on-line  [17], [12], [11] 
thanks  to  measurements.  However,  the  major  drawback  is  that  an 
LMS  controller  might  lead  to  a divergence  of  the  control  scheme  if, 
for  some  operating  conditions,  the  poles  of  the  IIR  become  unsta- 
ble. The  features  of  a LMS  controller  have  been  extensively  stud- 
ied by  Evesque  h Dowling  [12],[11],  and  here  it  was  necessary  to  in- 
troduce a parallel  algorithm  (based  on  the  Laguerre’s  method  [28]) 
to  prevent  a starting  divergence  due  to  the  controller.  Other  adap- 
tive schemes  already  developed  include  neural  networks  [6]  which 
is  a nonlinear  version  of  the  LMS-controlIer,  and  a minimisation 
scheme  based  on  the  downhill  simplex  algorithm  [25].  All  these 
schemes  provide  no  guarantee  that  the  controller  can  stabilise  the 
self-excited  combustion  process. 

An  efficient  way  to  prevent  any  divergence  of  the  adaptive  con- 
trol scheme  is  to  use  systematic  methods  for  designing  stable  adap- 
tive systems.  The  adaptive  controller,  called  STR  (Self-Tuning 
Regulator)  by  Annaswamy  et  al  [2],  is  designed  based  on  a Lya- 
punov stability  analysis  and  is  therefore  guaranteed  to  be  stable  for 


any  operating  conditions.  Furthermore,  the  STR  has  the  advan- 
tage of  avoiding  a system  identification  procedure  (which  is  one  of 
the  main  difficulties  in  implementing  a LMS  controller  [12]),  since 
it  uses  little  information  about  the  physical  process.  However,  so 
far,  only  a specific  simple  combustor  has  been  shown  to  have  the 
structure  required  for  the  design  of  a STR  [2].  Moreover,  the  STR 
cannot  accommodate  a time  delay  between  control  action  and  its 
detection. 

Our  purpose  in  this  paper  is  to: 

(i)  determine  the  general  features  of  a self-excited  combustion  sys- 
tem, rather  than  investigate  a particular  combustor  in  detail. 

(ii)  exploit  these  features  to  design  a novel  adaptive  controller 
that  is  guaranteed  to  stabilise  the  combustion  system,  the  major 
challenge  being  to  guarantee  stabilisation  in  the  presence  of  time 
delay  in  the  combustion  system. 

Step  (ii)  involves  first  the  choice  of  a low  order  fixed  controller 
structure  that  can  stabilise  the  system,  and  second  the  determi- 
nation of  adaptive  laws  for  the  controller  parameters  guaranteed 
to  converge  to  stabilising  values. 

Hence,  the  paper  is  divided  as  follows:  in  section  2,  a whole 
class  of  self-excited  combustion  systems  is  described,  and  its  fea- 
tures used  to  build  a controller  are  given.  Section  3 describes 
low  order  fixed  regulator  structures  that  can  control  a combustion 
system  with  or  without  time  delay,  while  section  4 deals  with  the 
design  of  a Self  Tuning  Regulator  (STR)  guaranteed  to  stabilise  a 
self-excited  combustion  process,  containing  or  not  some  time  de- 
lays. In  section  5,  the  performance  of  the  STR  is  illustrated  on  a 
simulation  based  on  a nonlinear  model  of  a premixed  ducted  flame 
developed  by  Dowling  [10]. 

2.  General  features  of  self-excited  combustion  systems 

Most  premixed  combustors  are  highly  resonant  systems  and 
may  develop  combustion  instabilities  for  some  operating  condi- 
tions. These  self-excited  oscillations  result  from  an  interaction 
between  unsteady  combustion  and  acoustic  waves:  unsteady 

combustion  generates  sound,  while  acoustic  waves  reflected  from 
the  boundaries  perturb  the  combustion  still  further.  Rather  than 
investigate  a particular  combustor  in  detail,  we  determine  the 
general  structure  of  this  interaction,  which  will  then  be  exploited 
to  design  fixed  and  adaptive  controllers  in  sections  3 and  4 
respectively. 


2.1  Open-loop  characteristics 

A wide  class  of  combustion  systems,  including  lean  premixed 
prevapourised  (LPP)  combustors  and  aeroengine  afterburners, 
can  be  modelled  as  a combustion  section  embedded  within  a 
network  of  pipes,  as  shown  in  figure  1.  We  will  investigate 
linear  low  frequency  perturbations  to  the  flow'  in  such  a pipework 
system.  The  flow  at  inlet  to  the  combustor  is  assumed  to  be 
isentropic,  and  the  frequencies  of  interest  are  low.  This  ensures 
that  the  combustion  zone  is  short  compared  with  the  wavelength. 
Moreover,  since  only  plane  waves  transport  acoustic  energy,  it  is 
sufficient  to  consider  one- dimensional  disturbances.  The  pressure 
and  velocity  upstream  the  flame  can  therefore  be  written  as  a 
linear  combination  of  the  waves  g and  /,  and  dowmstream  the 
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Figure  1:  Open-loop  self-excited  combustion  process  with  actuation 


flame  as  a linear  combination  of  the  waves  h and  j: 
in  — xu  < x < 0, 

P(M)  = P1+/(t-=^r)+g(t+ir^ir) 

u(xyt)  = u\  -f  r4r~  [/  (t  - z j-  ) -g(t+  - ^1  (l) 

P\C\  L \ Cl+UiJ  \ C\-U\J\ 

in  0 < x < xd) 

p{x,t)  = p2+h(t-^^)+j 

\ C2  + U2/  \ C2  — U2J 

u(x,t)  = U2  + Z7—  U (t-_  X_  j —j  \t+_  X-_  ) (2) 

P2c2  L \ C2+tX2/  \ C2-U2/J 

p,u,p  and  c are  respectively  the  pressure,  velocity,  density  and 
speed  of  sound.  The  suffices  1 and  2 denote  flow  quantities  just 
upstream  and  downstream  of  the  combustion  zone,  and  the  over- 
bar indicates  a mean  value.  The  distances  xu  and  Xd  are  the 
lengths  upstream  and  downstream  of  the  combustion  zone. 

The  boundary  conditions  of  the  combustor  are  characterized 
by  upstream  and  downstream  pressure  reflection  coefficients  Ru 
and  Rd  respectively.  Since  this  boundary  condition  neglects  the 
conversion  of  combustion-generated  entropy  waves  into  sound  at 
any  downstream  nozzle,  we  expect  it  to  be  a good  approximation 
when  the  time  taken  for  entropy  waves  to  convect  through  the 
straight  duct,  uzfxd,  exceeds  their  diffusion  time.  We  consider 
linear  disturbances  with  time  dependence  est , and  assume  that  in 
the  half  plane  Real(s)  > 0: 

\Ru(s)\  < 1, 

\Rd(s)\  < 1.  (3) 

Physically,  this  means  that  the  amplitude  of  the  reflected  wave 
is  less  than  the  incoming  wave  and  may  include  some  time  delay. 
Simple  duct  terminations  like  open  and  choked  ends  trivially  sat- 
isfy this  condition,  provided  appropriate  energy  loss  mechanisms 
are  included.  In  appendix  A,  we  show"  that  condition  (3)  is  also 
satisfied  by  reflection  from  general  pipework  configurations  with 
negligible  mean  flow.  The  reflected  waves  / and  j are  easily  ob- 
tained from  g and  h using  the  combustor  boundary  conditions* : 

/(£)  = Rn(s)[g(t  - ru)\  at  x = -xu 

j(t)  — «d(s)[/i(£  - rd)]  at  x = xdy  (4) 

where  ru  = 2xu/c\(l  — M\ 2)  and  rd  = 2xd/c2{l  — Ms2)  are  respec- 
tively the  upstream  and  downstream  propagation  time  delays,  and 
M is  the  mean  flow  Mach  number. 

*F(s}[  ] is  an  operator  of  the  variable  s = For  instance,  :r(i)  = 

;^y  [*/(*)]  means  that  ( t ) -f  x(t)  = y(t). 


The  equations  of  conservation  of  mass,  momentum  and  energy 
across  the  short  flame  zone  at  x = 0 can  be  written  in  a form  that 
is  independent  of  downstream  density  and  temperature  [9]: 


P2  — Pi  + PlUl(U2  — Ul)  — 0 

— f—{P2U2  ~pvn)  + -piui(ul  - u\)  ~ (5) 

7 — 1 l A 

Q is  the  instantaneous  rate  of  heat  release,  A is  the  combustor 
cross-sectional  area  and  7 is  the  ratio  of  specific  heat  capacities. 
Substitution  from  (1),  (2)  and  (4)  into  (5),  making  use  of  the  isen- 
tropic  condition  pi/p'l,  and  linearising  in  the  flow  perturbations 
give  the  time  evolution  of  the  outgoing  waves  g and  h generated 
by  the  unsteady  heat  release  Q(t): 

( Xu  X\2  \ ( g(t)  YnRu  YrzRd  \ ( [<?(*- t„)]  \ 

V *21  A'22  j V MO  )~\  Y2iRu  Y22Rd  ) V [h(t-rd)]  J 

+ ^ QUb-Q  ^ (6) 

where  Xij  and  Ytj  are  constant  coefficients  depending  on  the 
mean  flow  only,  and  given  in  appendix  B. 

After  taking  the  Laplace  transform1  of  the  system  (6)  and  us- 
ing (1)  and  the  boundary  condition  (4),  one  obtains  the  transfer 
function 


^ _ «,(s)  (RdYl2e~srd_Xi2XRue-STu_  l} 

G{s)  = qw  “ ( ) 

where 

( XU-  RuYi  X12  - RdYne~ST \ 

\ X21  - RuY2ie~ST*  X22  - RdY22G~~ST<i  J ' 

G(s)  describes  the  generation  of  unsteady  velocity  ui(t)  at  the 
flame,  due  to  the  unsteady  heat  release  Q(t), 

Since  the  self-excited  oscillation  results  from  a coupling  between 
unsteady  heat  release  and  acoustic  waves,  the  forcing  of  unsteady 
heat  release  due  to  incoming  flow  disturbances  at  the  flame  must 
be  also  described.  In  many  applications,  the  combustion  responds 
most  strongly  to  velocity  fluctuations.  This  is  because  in  acoustic 
waves  the  fractional  change  in  flow  velocity  is  order  M~l  larger 
than  the  fractional  change  in  pressure,  ie  a large  factor  at  the 
low  Mach  numbers  at  which  combustion  can  be  sustained.  This 
dependence  on  flow  velocity  can  either  be  seen  directly  through 
its  influence  on  flame  kinematics  and  shape  [10], [13],  or  indirectly 
through  its  influence  on  fuel-air  ratio  and  hence  on  the  rate  of 
combustion  in  LPP  systems  [29].  A transfer  function 

H(s)  = (8) 

Ul(s) 

is  introduced  to  describe  this  combustion  response.  In  many  cir- 
cumstances, H( s)  will  include  substantial  time  delays.  Models 
for  the  flame  transfer  function  H(s)  have  been  published  in  the 
literature  for  different  combustors.  However,  we  do  not  want  to  re- 
strict our  controller  design  to  any  particular  combustion  system  or 
model.  Instead  we  will  make  general,  non-restrictive  observations 
about  the  structure  of  the  transfer  function  H(s): 

• (I)  The  flame  is  stable  when  there  is  no  driving  velocity  u, 
which  means  that  the  poles  of  H(s)  are  ‘stable’  (ie,  are  in 
the  half  plane  Real(s)  < 0 and  so  lead  to  eigenmodes  with 
negative  growth  rate). 

* (II)  The  flame  response  has  a limited  bandwith,  therefore 
H - 4 0 when  s — > 00. 

These  assumptions  about  H(s)  fit  many  flame  models  given 
in  the  literature,  including  premixed  flames  [10]  and  LPP  sys- 
tems [14].  The  eigenfrequencies  can  be  determined  by  combining 
equations  (7)  and  (8).  They  satisfy 

tThe  same  notation  is  used  for  a temporal  signal  and  its  Laplace 
transform,  for  instance  u(t)  and  u(s). 
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1 -G(s)//($)  =0.  (9) 

When  a combustor  is  unstable,  there  are  roots  of  equation  (9) 
with  Real{s)>  0:  linear  perturbations  grow  exponentially  in  time. 
We  will  design  a feedback  controller  to  stabilise  such  a system. 

2.2  Actuated  combustion  system 


In  order  to  apply  active  control  to  our  self-excited  combus- 
tion system,  an  actuator  is  used  to  inject  a perturbation  and 
hence  break  the  damaging  coupling  between  unsteady  combustion 
and  acoustic  waves.  The  two  most  commonly  used  active  control 
inputs  are  loudspeaker  forcing  and  fuel-forcing  ( [18]  [26][4]).  We 
will  concentrate  on  fuel  forcing  which  is  the  most  relevant  for 
practical  applications.  An  actuator  is  driven  to  provide  extra 
fuel  (and  sometimes  air)  which  in  turn  produces  additional  heat 
release.  In  order  to  describe  the  impact  of  this  input  on  the 
combustion  system  characteristics,  we  study  the  relationship 
between  Vc,  a voltage  sent  to  the  actuator,  and  Pre/ , the 
fluctuating  pressure  measured  at  a location  xref  (see  figure  1). 
That  is,  our  goal  is  to  characterize  the  transfer  function 


If  the  unsteady  pressure  Pref  is  measured  upstream  the  flame 
(Xref  < 0),  then  Pref  is  a linear  combination  of  the  upstream 
waves  / and  g.  Using  (1)  and  the  boundary  condition  (4)  at 
x = — one  easily  obtains: 


-2  <(*u+*re/) 

P«/(a)  _ _ 1 + Rue 

zz  D\C\  

ui  ( s ) Rue~STv  — 1 


ec'i"ci 


(16) 


Therefore,  it  follows  from  (15)  and  (16)  that  the  open-loop  trans- 
fer function  of  the  actuated  system  given  in  (10)  can  be  written 
in  the  form: 


W(s)  = = W0(s)e-‘T““  , (17) 

Vc(s) 

where 


r tat-Tdet  + Ta  (18) 

is  the  total  time  delay  in  the  actuated  system,  ra  is  the  time  delay 
due  to  the  actuation,  = — a:re//(c i — u\)  is  the  detection  time 
delay  due  to  the  pressure  measurement  location,  and 


Vc(«) 


(10) 


which  represents  the  actuated  open-loop  combustion  process.  We 
derive  this  transfer  function  in  the  following. 

We  assume  that  the  fuel  injection  is  arranged  so  that  the  ex- 
ternal voltage  Vc  results  in  the  additional  heat  release  Qc  through 
the  following  transfer  function 


9M=Wa(s)e~^,  (11) 

Vc(s) 

where  Wa(s)  represents  the  actuator  dynamics.  Typically  the 
actuator  will  be  a valve  with  the  characteristics  of  a mass-spring- 
damper  system,  whose  dynamics  are  described  by  the  transfer 
function  Wa(s).  If  the  fuel-air  mixture  is  injected  directly  into  the 
combustion  zone,  the  combustion  response  will  be  instantaneous 
(ra  = 0).  However,  if  only  fuel  is  added,  there  will  be  a small  mix- 
ing time  delay  before  it  is  burnt  ( ra  > 0).  Often  it  is  hazardous  to 
inject  fuel  directly  into  the  flame.  If  the  additional  fuel  is  intro- 
duced some  distance  upstream  of  the  combustion  zone,  there  will 
be  a convection  time  delay  ra  between  injection  and  combustion. 
In  a LPP  system,  it  is  convenient  to  modulate  the  main  fuel  sup- 
plied in  the  premix  ducts,  in  which  case  ra  may  be  a significant 
proportion  of  the  period  of  the  self-excited  oscillations.  Notice 
that  rQ  is  independent  of  the  flame  radial  position,  which  means 
that  we  assume  the  same  time  delay  between  all  fuel  injection  and 
its  combustion.  When  the  combustion  zone  is  short  this  is  trivially 
satisfied.  If  the  combustion  zone  is  extensive,  it  may  be  necessary 
just  to  inject  fuel  in  a localised  region  to  meet  this  constraint. 

There  will  be  additional  unsteady  heat  release  driven  by  the 
flow  fluctuations.  We  will  denote  this  naturally  occuring  rate  of 
heat  release  by  Qn.  It  is  related  to  the  velocity  fluctuations  by 
the  flame  model  in  (8): 


H{s)  = (12) 

«l(«) 

For  linear  fluctuations,  we  can  superimpose  the  fluctuating  heat 
release  due  to  external  actuation  Qc  and  the  naturally  occuring 
heat  release  Qn  to  give  the  total  fluctuating  rate  of  heat  release: 


Q(5)  = QC(5)  + Qn(«).  (13) 

The  acoustic  waves  generated  by  Q{s)  are  described  by  (7),  ie 


G{b)  = 


m(s) 
Q(s)  ' 


From  equations  (11)- (13),  one  obtains  that: 


ui(s)  _ G(s)wa(s)e-ST» 
Vc(s)~  1 -G(s)H(s) 


(14) 


(15) 


Wo(s)  = Pi  Cl 


~2<(xu-t.g  f) 

G(s)Vya(s)  1 + ft„e 


1 -G(s)H(s)  - 1 

with  G(s)  and  H(s)  given  in  equations  (7)  and  (8). 


(19) 


If  Pref  is  measured  downstream  the  flame  ( xref  > 0),  then  Pref 
is  a linear  combination  of  the  downstream  waves  h and  j.  A sim- 
ple calculation,  using  the  wave  structure  in  equations  (1)  and  (2) 
and  the  continuity  condition  across  the  combustion  zone  in  equa- 
tion (5)a,  shows  that  the  transfer  function  Pref{s)jVc{s)  again 
has  the  form  given  in  (17),  provided  that  the  right-hand  side  of 
equation  (19)  is  multiplied  by 


D -Xn+RuYne-*r*  l+Rde 

Pdu(s)^~ 


2c  Xdr*r'l 

C2(  1— m|) 


Xn-RdYue-”0  *“**;;«/ 

1 + Rue 


U2-W2  ^ (20) 


2.3  General  structural  properties  of  the  open-loop  system 
useful  for  control  design 


For  the  sake  of  clarity,  only  the  case  of  an  upstream  pres- 
sure measurement  is  considered  here,  but  similar  results  can 
be  derived  for  a downstream  pressure  measurement,  using  the 
expression  of  Pdu(s)  given  in  (20).  The  wave  description  of  linear 
perturbations  in  sections  2.1  and  2.2  gives  the  open-loop  transfer 
function  W(s)  in  a particularly  convenient  form: 


^(5)  =e“srtotvy0(5)>  (21) 

ie  the  product  of  a pure  time  delay  and  the  transfer  function 
Wo(s)  defined  in  (19).  We  will  now  show  that  Wo(s)  has  some 
general  structural  properties  that  are  useful  for  the  control  de- 
sign. First  Wo  (5)  given  in  equation  (19)  is  expanded  into  a ratio- 
nal form.  This  is  achieved  by  applying  the  Fade  approximation 
technique  [3]  for  each  exponential  term  of  Wo(s).  This  technique 
has  been  widely  used  in  handling  systems  with  time  delays.  We 
will  use  the  notation  to  denote  the  (L,  M)th  order  Fade 

approximant  of  a function  f(s ),  which  is  a rational  function  P(s) 
whose  numerator  has  order  L and  denominator  order  M.  The  ra- 
tional function  P(s)  is  chosen  such  that  the  first  L + M + 1 terms 
in  the  power  series  of  P(s)  will  match  those  of  /(s),  ie 


m 


-}  + 0(sL+M+1). 

Mh(s) 


(22) 


Lemma  1:  The  zeros  of  Wo  are  stable  (ie  are  in  Real{s)  < 0): 
Wo(s)  is  said  to  be  ‘minimum  phase’ 
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Note  that  a rational  approximation  of  the  multiplying  time  delay 
e~srtot  in  equation  (21)  introduces  unstable  zeros  into  the  open- 
loop  transfer  function  W(s)  = Wo(s)e~sn°* : it  is  Wo(s)  that  is 
minimum  phase,  and  not  the  overall  transfer  function  W(s).  The 
time  delays  in  the  expression  of  Wo  (a)  given  in  (19)  arc  approxi- 
mated by  Pade  expansions  to  give  a rational  approximation  to  the 
transfer  function  IVo(s),  and  the  stability  of  the  zeros  of  Wo(s)  is 
discussed.  First',  it  is  evident  from  equation  (19)  that  the  poles 
of  the  flame  transfer  function  H(s)  become  zeros  of  Wo(s).  As 
noted  in  assumption  (I),  the  condition  that  the  flame  is  stable 
when  there  is  no  driving  velocity  ui  ensures  that  these  zeros  are 
all  in  Real(s)  < 0.  Therefore  H(s)  does  not  introduce  unstable 
zeros  for  Wq.  Secondly,  G(s)  appears  at  the  numerator  and  the 
denominator  of  Wo(s),  therefore  its  poles  do  not  influence  the  ze- 
ros of  Wb(s).  Thirdly,  we  assume  that  W0($)  has  no  unstable 
zeros.  Fourthly,  the  numerator  of  Wo  includes  G(s)  defined  in  (7) 
and  is  the  product  of  terms  of  the  form  K\(s)  4-  K2(s)e~3T , where 
r is  a time  delay,  and  | (s?)/ («) | > 1 for  Real(s)  > 0 (because 
the  reflection  coefficients  Ru(s)  and  Rd(s)  have  modulus  strictly 
smaller  than  1 for  Real(s)  > 0,  and  also  because  X12  > | V12 1) • 
These  factors  have  no  zeros  in  Real(s)  > 0,  because  at  a zero 


Ki  (f) 

K2(s)' 


(23) 


Equation  (23)  cannot  be  satisfied  in  the  half-plane  Real(s)  > 0, 
because  there  |e“5r|  < 1,  while  \Kx{s)/K2{s)\  > 1-  We  need  to 
check  that  this  remains  true  after  a suitable  Pad£  expansion  of 
e~ST . After  the  Pade  approximation  is  made,  one  has  to  solve 
equations  of  the  type 


L]  = -Ki  (s) 
MJ  K2(s)  * 


(24) 


in  order  to  find  the  remaining  zeros  of  Wq.  Baker  & Graves- 
Morris  [3]  introduce  the  concept  of  A-acceptability  for  rational 
functions:  a rational  function  R(z)  is  A-acceptable  if  |i7(,z)|  < 1 
for  Real(z)  < 0.  They  go  on  to  prove  that  the  Pade  approximant 
of  the  exponential  function  [L/M]ei  is  A-acceptable  provided  that 
M — L,  LTl  or  L 4-  2.  Therefore,  with  such  a choice  of  L and 
M,  we  obtain  that 


| [L/M]e-Jr  | < 1 for  Real(-ST)  < 0.  (25) 

However,  \K\ {s)f  Kz{s)\  > 1 for  Real(s)  > 0,  so  equation  (24) 
has  no  roots  in  Real(s)  > 0.  On  Real(s)  = 0,  the  numerator 
and  denominator  of  the  Pade  approximant  [L/M]e_,r  are  complex 
conjugate  when  L = M , which  means  that  |[M/M]e_J>r  | = 1 in 
Real(s)  — 0.  Hence,  since  \K\  ($)/K2(s)\  > 1 for  Real(s)  = 0, 
equation  (24)  has  no  roots  on  Real($)  = 0,  when  L = M. 

Therefore,  it  has  been  proved  that  all  the  zeros  of  Wq  are  stable, 
provided  that  the  (L,  M)th  order  of  each  Pade  approximant  is 
chosen  so  that  M = L,  for  any  M. 

Lemma  2:  The  relative  degree  of  Wq  is  equal  to  the  rela- 

tive degree  of  the  actuator  transfer  function 

An  essential  feature  required  for  the  control  design  is  the  relative 
degree  n*  of  Wo,  denoted  n*(Wo),  which  is  the  degree  of  the 
denominator  of  Wq  minus  the  degree  of  its  numerator.  n*(Wo) 
is  the  sum  of  the  relative  degrees  of  its  various  factors.  Again 
working  from  the  definition  from  the  definition  of  Wo  in  (19): 


• n (H)  > 0 since  the  flame  has  a limited  bandwidth  response 
(assumption  (II)).  Therefore,  H(s)  does  not  affect  the  relative 
degree  of  Wq  . 

Finally,  it  appears  that  the  relative  degree  of  Wq  is  equal  to 
the  relative  degree  of  the  actuator  transfer  function  Wa,  when 
each  Pade  approximant  [L/M]  satisfies  L < M. 

Lemma  3:  The  high  frequency  gain  of  Wq  is  positive 

The  controller  design  requires  information  on  the  sign  of  the  high 

frequency  gain  &o  of  Wo  [21],  defined  as  follows: 


W0(s)  = kQ 


Zq(s) 
Ro(s) ’ 


(26) 


where  ko  is  a constant,  and  Zq(s)  and  Rq(s)  are  two  monic* 
polynomials.  To  find  sign{ho),  we  simply  need  to  find  equivalent 
expressions  at  large  s for  each  factor  of  Wo-  As  noted  in  assump- 
tion (II),  H -4  0 for  s — > 00.  Furthermore,  n*(G)  = 0,  so  that 


1 - G(s)H{$)  ~ 1 for  s -+  00.  (27) 

The  other  factors  of  Wo  are  terms  of  the  form  14*  R(s)e~3T,  where 
r is  a time  delay.  Make  a Pade  approximation:  e~~ST  = [L/M], 
with  L = M.  Hence  the  high  frequency  gain  of  [L/M]  is  (— 1)M. 
In  appendix  A is  shown  that  the  high  frequency  gain  h of  R(s) 
satisfies  |h|  < 1.  Therefore,  the  high  frequency  gain  1 4-  h(  — l)M 
of  the  term  1 4-  R(s)e~ST  has  the  same  sign  as  1.  Finally,  at  high 
frequencies,  the  gain  of  Wo  is  easily  found  to  be  positive  when 
each  Pade  approximant  [L/M]  satisfies  L = M. 

There  is  a straightforward  reason  why  Wo(s),  in  the  open-loop 
transfer  function  PTtf/Vc  — Wo(s)e_sr<0< , has  the  simple  proper- 
ties outlined  in  lemmas  1 and  2.  For  Real(s)  > 0,  the  amplitudes 
of  the  oscillations  do  not  decrease  with  time  and  the  boundary 
conditions  |i?«|,]Aa|  < 1 ensure  that  the  largest  contribution  to 
PTef  is  from  the  acoustic  wave  leaving  the  combustion  zone,  rather 
than  the  waves  subsequently  reflected  from  the  boundaries.  Under 
these  circumstances,  the  main  structure  of  Wq(s)  in  equations  (19) 
and  (20)  is  dominated  by  the  properties  of  Wrt(s),  the  other  mul- 
tiplying factors  do  not  introduce  unstable  zeros  nor  affect  the  rel- 
ative degree. 

It  is  interesting  to  note  that  this  argument  remains  true  if  the 
form  of  actuation  is  a loudspeaker,  provided  the  loudspeaker  is 
located  within  the  combustion  zone.  However,  this  situation  is 
more  complicated  when  the  loudspeaker  is  at  a general  axial  po- 
sition in  the  combustor.  Then,  since  the  combustion  zone  is  an 
active  component,  it  can  reflect  a wave  of  greater  amplitude  than 
the  incident  wave:  if  Rf  denotes  the  reflection  coefficient  at  the 
flame,  \R/\  > 1 is  possible  even  in  Real(s)  > 0 [27].  Under  these 
circumstances,  Pref/Vc  can  have  unstable  zeros  for  some  positions 
Xrej ■ Lemma  1 therefore  is  not  true  for  general  loudspeaker  po- 
sitions. This  is  consistent  with  the  observations  of  Annaswamy 
et  al  [1]  who  calculated  Pref/Vc  for  a particular  idealised  com- 
bustor with  loudspeaker  actuation.  They  found  that  no  simple 
relationship  could  be  derived  between  the  locations  of  the  sensor, 
actuator  and  flame,  and  the  zeros  stability.  For  instance,  for  the 
particular  case  of  sensor  and  actuator  collocated  at  the  flame,  their 
open-loop  plant  Pre//Vc  had  no  unstable  zeros.  However,  for  the 
particular  case  of  fuel  forcing,  the  open-loop  transfer  function  sat- 
isfies lemmas  1 and  2,  properties  that  greatly  help  in  the  control 
design. 


• 7t*  (C7(s))  = 0 when  the  (L,  M)th  order  Pade  approximant  for 
each  exponential  term  e~$Tu  and  e~STd  satisfies  L < M (we 
assume  that  Ru(s)  and  Rd(s)  have  a relative  degree  equal  to 
zero,  as  it  is  described  in  Appendix  A). 

-2j(gu  +«  rtf) 

• similarly,  ri*{  e_1alu~_~1I } — ) = 0 when  the  (L,  M)th  or- 

der Pade  approximants  for  each  exponential  term  e~STu  and 

— 2j(xu  / ) 

e satisfies  L < M, 


3.  Fixed  regulator  design 

We  have  shown  that  for  fuel  actuation  the  open-loop  transfer 
function  has  a simple  form:  it  can  be  written  as  the  product  of  a 
pure  time  delay  and  a transfer  function  Wo(s)  which  is  rational 
after  a Pade  expansion  is  made:  Pref(s)/Vc(s)  = e~STtot  Wo(s). 
We  will  begin  by  designing  a controller  for  the  particular  case 

la  monic  polynomial  denotes  a polynomial  whose  leading  coefficient 
is  unity 
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Ttot  = 0 and  will  then  extend  these  ideas  to  the  general  and  more 
practically  relevant  case  Ttot  £ 0. 


left  half  plane  Real(s)  < 0,  However,  since  (35)  is  not  true  in 
general,  a better  strategy  is  to  use  the  following  regulator: 


3.1  System  without  time  delay  (rtot  = 0) 

It  is  clear  from  equation  (18)  that  Ttot  — 0 requires  that 
the  control  fuel  is  injected  and  burnt  with  no  time  delay  (ra=0), 
and  that  the  reference  pressure  is  measured  in  the  combustion 
zone  ( Xref  = 0).  Our  open-loop  combustion  process  is  then 
described  by 


Kj(s)  = fcci±ii,  (36) 

5 +Pc 

where  pc  and  zc  are  some  positive  constants,  pc  > zc . K'2(s) 
corresponds  to  a phase-lead  compensator,  which  adds  phase,  ie 
damping,  in  a frequency  range  [zCy  pc].  Then  the  closed-loop  poles 
are  the  zeros  of 


Wo  (s) 


Pref(s)  __  . Zq  ('S) 
Vc(s)  - °Ro{s) 


(28) 


where  ko  is  a positive  constant,  and  Zo  and  Ro  are  two  monic 
polynomials.  Furthermore,  Zo  and  Ro  are  ‘coprime’  polynomials, 
which  means  that  they  have  no  common  factors.  From  lemma  1, 
we  also  know  that  Z0  is  a stable  polynomial  (ie  it  has  only  zeros  in 
Real(s)  < 0),  whereas  Rq(s)  has  unstable  zeros  since  our  system 
exhibits  self-excited  oscillations.  Finally,  if  Rq  has  degree  n,  Zo 
has  degree  n — n* , where 


1 < n*  = n*(Wa)  < 2.  (29) 

(n*  = n*  (Wa)  comes  from  lemma  2,  and  most  practical  actuators 
have  a relative  degree  of  1 or  2). 

To  this  open-loop  system,  we  will  apply  an  active  feedback 


~ = ~K(S).  (30) 

Rrcf 

The  aim  of  the  regulator  is  to  stabilize  the  system,  ie  to  make  all 
the  closed-loop  poles  stable.  Combining  (28)  and  (30)  shows  that 
these  poles  are  the  zeros  of 


Rct(s)  = Ro(s)  + K(s)k0Z0(s ),  (31) 


Rci{s)  = (s  + Pc)Rq(s)  + kokc(s  + Zc)Zo(s).  (37) 

For  a ‘large  kcy  and  an  adequate  choice  of  zc  and  pc,  the  n zeros  of 
Rci(s)  can  be  moved  towards  the  left  half  plane.  More  precisely,  a 
finite  value  kCmin  > 0 and  some  positive  constants  pc  and  zc  exist 
such  that 


IM  ^ ^cm^n 

sign(kc)  — sign(ko) 

y ^(zerosof  Rq)  — (zeros  of  Zq)  < pc  — zc  (38) 

is  a necessary  and  sufficient  condition  to  stabilize  a minimum 
phase  plant  of  relative  degree  2 with  the  regulator  /^(s).  Notice 
that  K-i{s)  is  also  guaranteed  to  stabilize  a minimum  phase  plant 
of  relative  degree  1. 

In  the  following,  the  compensator  Kz(s)  given  in  equation  (36) 
will  be  implemented  as  shown  in  figure  2.  The  feedback  transfer 
function  for  this  system  is  given  by 

k2Vc 

Vc  = ~ " k\P ref  (39) 

S + Zc 

ie 

Vo  _ kljs  + Zcl  (40) 

P I ref  S -f  Zr  + ^2 


The  regulator  transfer  function  K(s)  is  to  be  determined  using  It  is  clear  from  equations  (40)  and  (36)  that  kx  represents  the 

roots  locus  arguments  [8]:  gain  kc  and  hz  determines  the  phase  lag  pc  in  iC>(s). 


If  tl*(Wo)  = 1,  consider  the  transfer  function 

Ki  (s)  = kc,  (32) 

where  kc  is  a constant.  Then  the  closed-loop  poles  are  the  zeros 
of  the 

Rods ) = Ro(s)  + kokcZo(s).  (33) 

For  )&c|  ‘large’,  n — 1 zeros  of  Rci(s)  will  be  moved  towards  the 
n—'  1 stable  zeros  of  koZo(s).  Investigation  of  the  large  \s\  asymp- 
totic form  shows  that  the  nth  zero  of  Rci(s)  is  also  stabilised  if 
sign(kc)  = sigrt(ko).  Therefore,  a finite  value  kCmin  > 0 exists 
such  that 

\kc\  > sign(kc)  = sign(ko)  (34) 

is  a necessary  and  sufficient  condition  to  stabilize  our  minimum 
phase  plant  of  relative  degree  1. 

If  n*(Wo)  = 2,  and  the  regulator  K\(s)  is  used,  then  a large 
\kc\  will  guarantee  that  n — 2 zeros  of  Rd(s)  will  be  moved  towards 
the  stable  zeros  of  k0Zo(s ).  The  two  remaining  complex  conjugate 
roots  of  Rd(s)  — 0 will  be  moved  towards  the  stable  half  plane 
Real(s)  < 0 only  if 

sign(kc)  — sign(ko) 

^^(zeros  of  Ro)  ~^2(zeros  of  Zo)  < 0.  (35) 

As  explained  by  Dorf  &c  Bishop  [8],  equation  (35) b guarantees 
that  the  asymptote  centroid  of  the  root  locus  is  situated  in  the 


Figure  2:  Fixed  low-order  controller  structure  for  Ttot  = 0,  n*  <2 

However,  the  major  drawback  of  such  a fixed  regulator  Kz{s) 
is  that  a cautious  choice  of  pc  is  necessary  if  the  inequality  (38)c 
is  to  be  guaranteed  without  detailed  knowledge  of  the  plant.  This 
in  turn  can  mean  that  the  gain  kc  required  to  achieve  control  is 
large,  especially  to  ensure  stabilization  under  varying  operating 
conditions,  ie  under  uncertainties  in  the  unstable  frequencies  u?u  • 
A large  kc  means  a large  control  effort,  which  is  to  be  avoided. 
Therefore,  to  improve  the  response  of  our  regulator  Kz(s)  under 
varying  operating  conditions,  one  can  choose  a fixed  zc  > 0,  and 
make  the  other  control  parameters  kc  and  pc  adaptive.  This  is 
the  topic  of  section  4 . 

3.2  Combustion  system  with  known  time  delay  (rtot  ^ 0) 

Here  the  combustion  process  is  described  by 

W(s)  = = k0  = WoWe-"*- , (41) 

Vc($)  ^o(s) 

with  ko  a constant,  Zo  and  Ro  two  coprime  and  monic  polyno- 
mials, and  not  is  a known  time  delay.  Ro  has  degree  n,  and  Zo 
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has  degree  n — n* , with  n*  = 1 or  2 (lemma  2).  It  was  shown  in 
section  2 that  Z0  is  a stable  polynomial  (lemma  1)  and  that  ko  is 
positive  (lemma  3). 

In  section  3.1,  low  order  controllers  of  the  form  (32)  and  (36) 
were  used  to  control  the  combustion  process.  Here,  the  presence  of 
the  time  delay  rtot  makes  such  controllers  inadequate,  especially 
when  rtot  is  of  the  order  of  the  period  of  the  unstable  mode.  Con- 
trol of  systems  in  the  presence  of  time  delays  has  been  extensively 
studied  [30] [19] [15].  A popular  approach  to  accommodate  large 
delays  is  due  to  O.J.M.  Smith  [30].  The  Smith  Controller  (SC) 
attempts  to  estimate  future  outputs  Prf>f  of  the  system  using  a 
known  model,  and  provides  an  appropriate  stabilization  action. 
In  M an i tius  & Olbrot  [19],  the  SC  was  modified  to  control  sys- 
tems that  are  open-loop  unstable  by  using  finite-time  integrals  of 
inputs  Vc  to  estimate  future  outputs  (the  goal  of  this  modification 
was  to  avoid  unstable  pole-zero  cancellations).  In  Ichikawa  [15], 
a SC  controller  with  finite-time  integrals  was  derived  to  serve  as 
a pole-placement  controller.  In  [24]  and  [16],  adaptive  versions  of 
the  pole- placement  controller  have  been  developed  and  proved  to 
be  stable. 

Since  our  plant  is  open-loop  unstable,  we  want  to  implement 
a SC  using  the  finite-time  integral  suggested  by  Manitius  Sz  Ol- 
brot [19]  and  used  in  Ichikawa  [15],  and  which  is  given  by 

o 

Vsa(t)~  j A (a)Vc{t  + <T)da,  (42) 

— Ttot 

where  A(<r)  is  a weighting  function  (it  is  shown  in  section  4.2  how 
discrete  values  of  A(.)  can  be  found  by  adaptation). 

Briefly,  the  controller  used  in  [15]  and  shown  in  figure  3 has  the 
following  form: 

o 

VcW  = -|[^[V,c(i-r(<>e)1-|M[Fre/(t)]+  J A (a)Vc(t+a)d<r  (43) 

-rtot 

where  the  operators  C(s)}  D(s)  and  E(s)  are  chosen  such  that  the 
denominator  of  the  closed- loop  transfer  function  matches  a chosen 
stable  polynomial.  This  can  be  achieved  only  if  degree(E)  = n, 
degree{C)  = n — 2 and  degree(D)  — n — 1.  This  means  that  even 
for  a plant  of  low  relative  degree  (n*  = 1 or  2),  if  the  order  n is 
large,  then  the  controller  dynamics  C/E  and  D/E  will  be  of  high 
order.  This  differs  from  the  delay-free  case  described  in  section 
3.1,  for  which  a first  order  compensator  could  stabilize  the  plant 
of  relative  degree  2 or  1. 


Figure  3:  A nth  order  controller  structure  for  rtot  ^ 0,  given  in 
Ichikawa’s  paper  [15] 

However,  our  aim  is  to  keep  a low  order  compensator  to  control 
the  plant  as  it  was  done  in  section  3.1.  Therefore,  we  suggest 
the  controller  structure  given  in  figure  4.  We  begin  as  in  [15]  by 
choosing  A(cr),  the  weighting  function  in  the  Smith  Controller,  to 
have  the  form 

A(cr)  = (44) 

*=  1 

where  pi  are  chosen  to  be  the  n zeros  of  the  polynomial  Rq(s ), 
whereas  (at  this  stage)  the  a,  are  arbitrary  coefficients.  Then,  af- 
ter substitution  into  (42)  and  evaluation  of  the  integral,  we  obtain 


where 

(45) 

ni(5)  _ A Oi 

Rq(s)  iC\s~ 

n2(s)  = y ' eftr...  (46) 

Ro{s)  5 ~ Pi 

We  prove  in  the  following  that  the  low  order  controller  associated 
with  the  SC  (see  figure  4)  will  stabilize  our  minimum  phase  plant 
of  relative  degree  1 or  2 for  a small  rtot.  With  the  controller 
structure  described  in  figure  4,  the  closed-loop  transfer  function  is 
given  by 

where  the  closed-loop  poles  are  the  zeros  of 

, (47) 

Rd(s)  = A(s)  ^ B(s)e~sr^  y 

with 

(48) 

A(s)  = ( s 4 zc)(Rq(s)  - + k2R^(s) 

B(s)  = (s  4 zc)  (n2($)  4 k\koZo(s)) . 

(49) 

(50) 

Figure  4:  A fixed  low-order  controller  structure  for  rtot  # 0,  n*  <2 


If  rtot  = 0,  then  equation  (46)  gives  that  r*i(s)  = 712(5).  There- 
fore, 

Rci{s)  = {s  + zc  + k2)Ro(s)  -I-  k\ko(s  4-  zc)Zq(s) 

= M(s)  say.  (51) 

That  is,  the  closed-loop  polynomial  coincides  with  that  in  equa- 
tion (37),  which  was  for  the  delay-free  case.  As  shown  in  section 
3.1,  Ao(s)  is  stable  for  certain  choices  of  k\  and  &2. 

If  Tt0t  0,  we  show  that  the  controller  structure  given  in  fig- 
ure 4 still  guarantees  stability.  Since  the  polynomial  712(5)  is  of 
degree  n — 1,  its  coefficients  can  be  chosen  such  that 

7*2(5)  = — kikoZo(s))  (52) 

which  implies,  from  equation  (50),  that  B(s)  = 0 and  therefore 
Rcl(s)  = .4(s). 

Furthermore,  it  is  clear  from  equation  (46)  that  the  polynomials 
m and  «2  of  the  SC  are  linked.  More  precisely,  the  choice  of  7*2 
in  (52)  imposes  restrictions  on  ni : it  means  that  the  coefficients  aM 
and  hence  ni,  are  proportional  to  koki . When  rtot  = 0,  n 2 = 7ii, 
we  will  emphasise  this  scaling  by  writing: 

7*2(5)  = 7*1  (s)  4-  rtotkokin^s),  (53) 

where  ro3  is  a polynomial  of  degree  n — 1,  with  finite  coefficients. 
Therefore,  using  (52)  and  (53)  in  (48),  the  closed  loop  poles  are 
the  roots  of 

Rcl{s)  = (s-f  Zc-hk2)Ro{s)+kiko{s  + zc)  [Zo(tf) -n<>trc3(s)]  (54) 
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For  ‘small1  rtot,  ie  for  rt0twa  <0(1)  where  uja  is  the  highest 
frequency  among  the  zeros  of  Z0,  the  zeros  of  T(s)  = Zq(s)  - 
Ttotns(s)  are  close  to  the  zeros  of  Zo(s)i  and  hence  are  stable. 
Therefore,  for  |Aq|  large,  n — 1 zeros  of  Rci{s)  are  stabilised.  The 
2 remaining  zeros  of  Rci{s)  are  obtained  at  large  $,  (k^l(s  + zc-f 
k2)Ro{s)  is  not  negligible  compared  to  (s  + zc)(Z0(s)  - Ttotm(s)) 
when  s and  this  is  where  we  find  the  2 remaining  zeros 

of  Rd(s)).  After  division  by  sn~ 1 of  the  3 highest  coefficients  of 
Rci(s ),  we  obtain 

s2  + (k‘2  — kikoTtotC)s  + (kiko  ~ kikoTtotCzc)  (55) 

where  C is  the  highest  coefficient  in  n3(s).  The  2 remaining  roots 
of  Rci(s ) are  stable  if  the  polynomial  (55)  is  stable,  ie  if 

k2  — kikortotC  > 0 

k\ ko(l  — rtotCzc)  > 0.  (56) 

In  the  range  of  values  of  not  considered  (ie  TtatUa  <0(1)),  one 
easily  checks  that  C~+ 1 and  that  (56)  is  satisfied  for  some  k\  >0 
and  k2  > 0.  Therefore,  we  proved  that  for  ‘small  all  the 

closed-loop  poles,  ie  the  roots  of  iM(s),  are  stabilised  for  some 
ki  > 0 and  k2  > 0.  In  other  words,  a plant  of  relative  degree  nm  <2 
with  a time  delay  rtot.  not  too  large  is  stabilized  by  the  controller 
structure  given  in  figure  4,  and  the  controller  equation  is  given  by: 

o 

Vc(t)  = -kiPref(t) — [V^(t)]  + f A(a)Vc(t  + a)  da.  (57) 

s + zc  J 

~Ttot 


•* — * — * < — * * — * 

tmtful  (tbuiincd  contrul  control  control  control 
in  simulation 


Figure  5:  Comparison  of  simulation  results  for  an  infinite  order  plant 
(model  described  in  section  5,  ufu  = 370  rad/s)  with  theoretical  results 
for  a second  order  plant  (appendix  C).  (a)  k\  obtained  in  simulation  and 
sign(k\)  required  for  stability  for  a second  order  plant,  (b)  simulation: 
control  obtained  periodically  up  to  rtot  ~ 3 cycles  of  oscillations. 


In  practice,  the  constraint  on  the  size  of  n0t  is  not  so  strong: 
simulation  results  for  a nonlinear  model  of  an  infinite  order  plant 
descibed  in  section  5,  show  that  control  is  obtained  with  k\  > 0 
and  k2  > 0 for  rtot  < where  u ;u  is  the  main  unstable  mode. 
Furthermore,  for  higher  values  of  Ttot , up  to  n0t  ~ 3 cycles  of  oscil- 
lations, we  observe  some  periodic  stability  bands  according  to  the 
values  of  ujurt0t  (figure  5b):  a ‘stability  band’  corresponds  to  val- 
ues of  utuTtot  for  which  control  is  obtained  after  a finite  time  called 
settling  time.  Between  two  consecutive  ‘stability  bands’,  there  are 
a few  values  of  ujun0t  for  which  control  is  not  obtained  (the  set- 
tling time  is  then  infinite).  It  was  also  observed  that  the  sign 
of  the  first  order  compensator  gain  k\  achieving  control  changes 
between  two  consecutive  stability  bands  (see  figure  5a).  These 
observations  on  sign(k\)  and  on  the  stability  bands  pattern  can 
be  interpreted  as  follows:  for  frequencies  smaller  or  equal  to  cuu, 
the  open- loop  transfer  function  can  be  approximated  by  a second 
order  system 


*o 


(5  - <Tu)2  + <*>£ 


(58) 


where  k'0  = k0(*l  is  a positive  gain*. 

It  is  shown  in  appendix  C that  a plant  of  order  2,  whose  transfer 
function  is  given  by  (58),  is  stable  for  any  delay  rtot  which  satisfies 


4.  Adaptive  regulator  design 

4.1  System  without  time  delay  (rtot  = 0) 


It  was  shown  in  section  3.1  that  a first  order  regulator 
jftT2(*s)  — Ms  + zc)/(s  + pc)  will  stabilize  our  combustion  process 
which  is  minimum  phase  and  of  relative  degree  less  or  equal  to 
2.  The  adaptive  version  of  this  regulator,  called  Self- Tuning 
Regulator  (STR),  is  given  in  figure  6. 


Figure  6:  Low  order  adaptive  controller  for  rtot  = 0,  n*  <2 


With  this  controller  structure,  the  closed-loop  transfer  function 
between  some  input  noise  i and  Pref  is 


ki  si n(o?u  rtot)  < 0 

1/tan (wurtot)  - ~ < 0.  (59) 

UJu 

Hence,  such  a plant  is  characterized  by  a pattern  of  periodic 
stability  bands  according  to  the  values  of  ujurtot,  as  shown  in  figure 
5a.  We  see  also  from  (59)  that  the  sign  of  ki  required  for  control 
satisfies 


sign(ki)  = —$ign(s\n{u>uTtot)),  (60) 

which  corresponds  to  observations  from  the  simulation  (see  figure 
5a).  In  a practical  combustor,  the  time  delay  not  would  usually 
not  exceed  3 cycles  of  oscillations,  therefore  our  controller  appears 
very  adequate  for  control  of  combustion  oscillations. 

§the  zeros  of  Rq  are  complex  conjugates,  hence  the  constant  coeffi- 
cient of  Rq  is  positive.  ZQ  is  monic  and  stable,  therefore  its  constant 
coefficient  is  also  positive. 


Prefjs) 


= Wct(s)  = 


kp(s  + Zc)Zq(s) 


<W  Rd(s) 

and  the  closed-loop  poles  axe  the  roots  of 


(61) 


fAi(s)  = (^  + zc  4-  k2)Ro(s)  + fciMs  + Zc)Z0{s),  (62) 

where  zc  > 0 is  fixed,  and  two  controller  parameters,  ki(t)  and 
M*)>  are  tuned.  It  is  clear  from  equations  (62)  and  (37)  that 
k\  represents  the  gain  kc  of  the  fixed  regulator  iMs),  while  k2 
is  used  to  tune  pc.  In  the  following,  vectors  are  denoted  in  bold 
characters  and  T denotes  the  transpose  of  a vector.  We  introdxice: 

* the  unknown  controller  parameter  vector  k (t)T  — 
[ — — fc2(t)],  and  the  error  parameter  vector  k = k-k", 
where  * denotes  a value  for  which  closed-loop  stability  is 
achieved. 
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♦ the  data  vector  d(t)T  = [Pre/(t),  where  V(t)  = 

j^[Vc(t)]  and  zc  is  a positive  constant.  In  an  experiment, 
d can  be  determined  at  each  t from  the  measurement  Pref 
and  known  Vc. 

We  now  need  to  find  an  updating  rule  for  k,  so  that  stabilisation 
can  be  accomplished  fo  any  parameters  in  Wo(s). 

Case  (i)  : ti*(Wq)  = 1.  This  corresponds  to  a case  for  which 
Narendra  &:  Annaswamy  [21]  have  developed  a STR.  However, 
we  will  repeat  the  main  points  of  their  argument  because  they 
provide  the  background  to  our  novel  STR  for  the  case  with  time 
delay  in  section  4.2.  When  n*  (Wo)  = 1,  the  closed-loop  transfer 
function  Wci  given  in  (61)  has  then  a relative  degree  equal  to  1 and 
is  ‘Strictly  Positive  Real’  (SPR),  which  is  the  essential  property 
required  to  develop  a global  stability  analysis  based  on  Lyapunov’s 
direct  method  [21].  Such  a method  aims  at  finding  adaptive  laws 
for  k which  are  guaranteed  to  stabilize  the  self-excited  combustion 
process.  Essentially,  the  strictly  positive  realness  of  Wci  means 
that  it  is  possible  to  find  a quadratic  positive  function  Vj,  that 
decays  in  time  when  k is  updated  correctly.  Such  a function  V\ 
is  referred  to  as  a ‘Lyapunov  function’,  and  can  be  viewed  as  an 
energy  function:  if  this  function  decreases,  it  implies  that  the 
system  is  stabilised.  To  summarize  Narendra  & Annaswamy’s 
results  [21],  when  Wo  has  a relative  degree  equal  to  1,  the  STR 
which  guarantees  the  stability  of  the  system  is  described  by: 

Vc(t)  = kT(t).d(t) 

k(t)  = -sign(k0)Prcf(t)d(t),  (63) 

where  sign(ko)  = +1  from  Lemma  3.  For  reference,  equa- 
tion (63)  comes  from  the  application  of  lemma  5.1  in  Narendra 
& Annaswamy  [21],  noting  that 

Pref(t)  = V^cK«)[kT(0-d(0]*  (64) 


Case  (ii)  : ra*(FFo)  = 2.  In  this  case  Wci  is  not  SPR,  but 
the  approach  of  Annaswamy  et  al  [2]  shows  that  modifying  the 
control  signal  Vc  as  indicated  in  figure  7 effectively  makes  the 
closed-loop  transfer  function  Wci  have  relative  degree  1 . Following 
their  approach,  we  write 

Vc(t)  = (s  + a)[kT(t).da(t)]  (65) 

where  da(£)  = j^[d(t)].  After  some  simple  algebra,  Vc  can  be 
written  as 

Vc(t)  = kT(t).d(t)  + kr(£).da(t).  (66) 

Now,  we  have 

Pre.f(t)  = VFC|(s)($  + a)[kT(t).da(t)] 

= Wm(5)[kT(t).da(t)],  (67) 

where  Wm(s)  = (s  4-  a)Wci(s)  is  the  new  closed-loop  transfer 
function  of  our  system,  which  has  relative  degree  2 and  is  SPR. 
Noting  that  equation  (67)  looks  similar  to  equation  (64),  lemma 
5.1  [21]  can  be  applied:  the  STR  which  will  guarantee  the  stability 
of  the  system  of  relative  degree  2 is  given  by: 

Vc  (t)  = kT(£).d(t)  + k(£)T.da(t) 

k(£)  = -sign(k0)Pr€f(t)dri(t),  (68) 

where  sign(ko)  = H-l  from  Lemma  3. 

An  important  remark  is  that  the  STR  described  for  the  case 
n*(Wo)  = 2 will  also  work  if  n*(Wo)  = 1,  which  means  that  even 
if  the  actuator  dynamics  are  not  very  well  known  in  practice,  the 
STR  design  given  in  equation  (68)  will  always  give  satisfactory 
results. 


Figure  7:  Modification  of  the  control  input  Vc 


4.2  Combustion  system  with  known  time  delay  ( rtot  ^ 0) 

The  controller  structure  given  in  figure  4 includes  fixed  controller 
parameters  kit  k% , a,  and  pi  which  need  to  be  chosen  based  on 
the  system  parameters.  Under  uncertainties  and  variations  in 
the  operating  conditions,  it  is  more  appropriate  to  adapt  those 
control  parameters.  Hence,  we  choose  an  adaptive  controller 
structure  as  shown  in  figure  8. 


Figure  8:  A low-order  adaptive  controller  for  Ttot  ^ 0 , n*  <2 

In  the  control  law  given  in  equation  (57),  the  finite-time  inte- 
gral due  to  Smith  Controller  and  given  in  equation  (42)  is  approx- 
imated as  follows: 

N 

Vsc(t)  = Ai(t)V«(t  - idt).  (69) 

t=l 

Similar  to  the  delay-free  case,  we  define  the  controller  parameters 
and  data  vectors  k and  d respectively: 

• k(t)T  = [— &i(t),  — &2(£),  Av(f), •••,  Ai (£)],  and  its  error  vector 
k = k-k\ 

• d (t)T  = \Pref(t),V(t),  Vc(t  - Ardt),...,  Vc(t  - dt)],  where 

m = TfciKcWj. 

Therefore,  as  the  time  delay  rtoe  is  increased,  N must  be  in- 
creased and  more  controller  parameters  are  required.  We  are  look- 
ing for  an  updating  rule  for  k. 

General  case  : n*(W0)  <2.  As  in  section  4.1,  the  closed 
loop  transfer  function  Wci0  given  in  (47)  is  made  to  effectively 
have  relative  degree  1 by  using  a modified  control  signal  Vc  (see 
figure  7).  Then  we  obtain: 

Prc/M  = VtW'X*  + [fcT(t).<U(01,  (*>) 

where  d a(f)  = ~-[d(t)].  Due  to  the  presence  of  the  time  delay 
e~STtot  ^ ^jie  LyapUnov  function  Vi  used  in  the  delay- free  case  will 
not  decay  in  time.  So  we  need  to  add  an  extra  positive  term 

o t j. 

in  Vb  in  the  form  of  a double  integral  j f ||k(£)||2d£di/,  to 

—rtoi  t+v 

account  for  the  time  delay  rtot.  Then  we  have  V\  < 0,  and  hence 
a system  whose  energy  is  decaying  in  time  (see  details  of  proof 
in  Appendix  D).  Therefore,  with  the  Lyapunov  function  given  in 
equation  (98),  the  stability  of  the  system  is  guaranteed  when  one 
uses  the  following  STR: 

Vc(t)  = kT(t).d(t)+k(t)TMt) 

k(£)  = -sign(k0)Pref(t)dtl(t -Ttot),  (71) 

where  sign(ko)  — +1. 
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Particular  case  : n*(VFo)  — 1.  Then  the  manipulation  on 
Vc  shown  in  figure  7 is  not  required,  and  therefore  the  following 
simpler  algorithm  can  be  implemented: 

Vc{t)  = kT(t).d{t) 

k(t)  = -sign(k0)PTef(t)d(t  - Ttot),  (72) 

where  sign(ko)  — +1. 

5.  Application  to  a premixed  ducted  flame  (simulation) 

5.1  Nonlinear  flame  model  used  in  the  simulation 

A model  for  nonlinear  oscillations  of  a ducted  flame  devel- 
oped by  Dowling  [10]  will  be  used  to  verify  the  controller 
performance.  The  theory  involves  extension  of  the  flame  model 
of  Fleifil  et  al[13]  to  include  a flame  holder  at  the  centre  of  the 
duct  and  nonlinear  effects.  Due  to  page  limitations,  the  reader  is 
asked  to  refer  to  [10]  for  details. 

5.2  Adaptive  regulator  design 


For  both  cases,  a convergence  coefficient  p is  added  in  the 
adaptive  law  for  each  controller  parameter.  The  controller 
parameter  vector  k is  initialized  to  zero.  Simulations  results  are 
given  in  the  next  section. 

5.3  Simulation  results 


System  without  time  delay  (rtot  — 0)-  Figure  9 shows 
the  time  evolution  of  the  pressure  measurement  Pref  and  the  cor- 
responding control  signal  Vc  for  varying  operating  conditions.  The 
control  is  switched  on  at  t = 0.15  s,  when  the  limit  cycles  are  al- 
ready established  (operating  conditions:  cp  = 0.7,  Mi  = 0.08). 
The  oscillations  tend  to  zero  within  1 s;  note  that  the  settling 
time  depends  on  the  convergence  coefficient  ft  used.  Then  from 
t = 1.95  s to  t = 2.1  s,  Mi  is  increased  to  0.095:  in  the  open-loop 
system  such  a change  shifts  the  frequency  of  the  unstable  mode 
from  58  to  63  Hz.  Nevertheless,  the  STR  is  able  to  maintain  con- 
trol. Changes  in  the  fuel  air  ratio  <p  have  also  been  successfully 
tested.  It  has  been  observed  that  the  STR  works  fine  even  if  xrej 
is  chosen  anywhere  upstream  the  flame,  which  is  interesting  for 
practical  applications  where  a pressure  measurement  directly  at 
the  flame  may  be  difficult  to  perform. 


The  flame  model  described  in  [10]  fits  in  the  class  of  com- 
bustion systems  given  in  section  2,  since: 

• the  upstream  reflection  coefficient  has  a modulus  strictly 
smaller  than  1 (choked  end:  Ru  = (1  ~ Mi) /(l  + Mi))  , in 
agreement  with  equation  3.  Notice  that  the  downstream  re- 
flection coefficient  has  modulus  just  equal  to  1 (ideal  open 
end:  R j — — l).  However,  with  care  the  condition  |Rd(«s)|  < 1 
in  equation  3 can  be  relaxed  to  |Rd(s)|  < 1 provided  an 
upstream  pressure  measurement  is  made.  This  is  because 
in  the  transfer  function  G(s)  given  in  (7),  IA12I  > IW2I  and 
IX22I  > \Y22l 

• after  linearisation  for  small  perturbations,  the  flame  transfer 
function  is  given  by  (see  [10]): 


H(s)  = 


Qn($) 

uG(s) 


2r}AHe-3T 2 

sTfcl  (6  + a) 


, q~-  . b CL  , ST  t\ 

a— be  8 * -I (1— e f) 

ST f 


(73) 


Therefore,  H 0 as  s -*  00,  and  H(s)  has  no  poles:  hence 
assumptions  (I)  and  (II)  are  satisfied. 


Hence,  section  2 proves  that,  once  the  system  is  approximated 
as  finite  dimensional,  the  relative  degree  of  the  open-loop  transfer 
function  kFo  is  equal  to  the  relative  degree  of  the  actuator  transfer 
function.  As  in  Evesque  et  al  [11],  the  actuator  chosen  is  a fuel 
injection  system  modelled  as  follows: 


= (74) 

V^(s)  s2/ta|  + 2 cs/ljc  + 1 

where  ujc  and  c are  respectively  the  resonance  frequency  and 
damping  of  the  fuel  injector.  From  equation  (74),  we  deduce  that 
the  relative  degree  of  Wo  is  2.  Section  2 shows  also  that  the  ‘high 
frequency  gain’  ko  is  positive  and  that  the  zeros  of  Wq  axe  stable. 

Two  cases  are  studied  in  the  simulation: 


Figure  9:  STR  for  = 0,  Ml  is  varied  linearly  from  0.08  at  t = 1.95 
s to  0.095  at  t — 2.1  s.  ...:  controller  OFF,  — : controller  ON  after 
t = 0.15  s. 


System  with  time  delay  ( rtot  7^  0).  This  is  the  most  in- 
teresting case  for  practical  applications.  As  already  illustrated  in 
section  3 (see  figure  5),  control  is  obtained  for  values  of  rtot  much 
larger  than  those  predicted  by  the  theory:  control  is  achieved  pe- 
riodically up  to  rtot  ~ 3 cycles  of  oscillations  (ie  up  to  uuTtot  ~ 67t 
where  is  the  unstable  mode).  The  settling  time  varies  peri- 
odically: control  is  easier  in  the  middle  of  a stability  band,  and 
becomes  more  difficult  on  the  boundaries  of  a stability  bands.  As 
in  the  delay  free  case  {rtot  — 0),  the  actual  value  of  the  settling 
time  depends  on  the  convergence  coefficient  ft  chosen.  A typical 
time  evolution  of  the  pressure  oscillations  while  control  is  applied 
is  shown  in  figure  10.  There,  in  order  to  test  the  adaptability  of 
our  STR  under  varying  operating  conditions,  the  mean  upstream 
Mach  number  Mi  is  increased  from  0.08  to  0.095  between  t = 1.95 
s and  £ = 2.1  s:  the  STR  successfully  maintains  Pref  at  zero. 

6.  Conclusions 


• a system  without  time  delay  (rtot  = 0).  This  is  achieved  by 
making  the  pressure  measurement  at  the  flame  (ie  xre/  = 0) 
and  by  setting  the  time  delay  ra  to  zero  in  equation  (74). 
Then  the  STR  is  implemented  as  indicated  in  equation  (68), 
with  sigTi(ko)  = +1. 

• a system  with  time  delay  ( rtot  ^ 0).  The  pressure  mea- 
surement is  chosen  for  instance  at  xrej  = — xu/2  and  the 
delay  ra  can  be  varied.  Then  the  adequate  STR  to  control 
the  self-excited  oscillations  is  given  in  equation  (71),  with 
sign(ko)  = +1. 


A general  combustion  system,  susceptible  to  combustion  instabil- 
ities, and  satisfying  the  following  assumptions,  was  considered: 

L the  pressure  reflection  coefficients  at  the  boundaries  of 
the  combustor  have  a modulus  strictly  smaller  than  1 in 
Real(s)  > 0. 

2.  the  flame  response  has  a limited  bandwidth:  at  high  frequen- 
cies, the  flame  does  not  follow  incoming  velocity  fluctuations. 

3.  the  flame  is  stable  in  itself,  it  becomes  unstable  only  due  to 
the  interaction  with  the  acoustic  waves  in  the  combustor. 
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Figure  10:  STR  for  rtot  = 23  ms,  xref  = -~xu/2,  M i is  varied  linearly 
from  0.08  at  t = 1.95  s to  0.095  at  t = 2.1  s.  . . . : controller  OFF,  — 
: controller  ON  after  t ~ 0.15  s. 


4.  the  actuator  used  for  control  is  a fuel  injector  which  produces 
a heat  release  rate  Qc  related  to  the  voltage  Vc  driving  the 
actuator  simply  by  a time  delay  and  a first  or  second  order 
differential  equation. 

We  demonstrated  that  such  a combustion  system  is  essentially 
represented  by  an  open-loop  transfer  function 


W(s) 


P ref (s) 

Vc{s) 


= ko 


Zq (£)  _. 
Ro  (s)€ 


where  Pre/  is  a pressure  measurement  in  the  combustor  and  Vc 
is  the  voltage  driving  the  fuel  injector  used  as  control  actuation. 
Some  general  properties  have  been  derived: 

(i)  the  zeros  of  Zo(s ) are  all  stable  (ie  are  in  Real(s)  < 0) 

(ii)  Once  Zo/Ro  has  been  made  rational  using  a Pade  expansion, 
the  relative  degree  of  Z0(s)/Ro(s ) is  equal  to  the  relative  degree 
of  the  actuator  transfer  function,  ie  2 for  a fuel  injector. 

(iii)  ko  is  a positive  gain. 

These  properties  have  been  exploited  to  design  an  adaptive  con- 
troll er  guaranteed  to  stabilise  the  self-excited  combustion  system. 
In  particular,  for  the  case  of  a combustion  system  with  time  delay 
(rtot  7^  0),  which  is  the  most  realistic  case,  the  design  is  completely 
novel:  it  involves  a first  order  compensator  combined  with  a Smith 
Controller.  The  adaptive  laws  for  the  controller  parameters  are 
derived  based  on  a Lyapunov  stability  analysis.  The  adaptive  con- 
troller has  been  tested  successfully  on  a simulation  of  a premixed 
ducted  flame.  An  experimental  verification  is  planned  over  the 
next  months. 
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APPENDIX 

Appendix  A.  Consider  linear  disturbances  with  time  depen- 
dence est  in  the  pipework  system  upstream  the  combustion  zone, 
as  shown  in  figure  11.  The  Mach  number  considered  is  so  low  that 
we  assume  there  is  no  mean  flow  in  the  pipework  system.  The  pres- 
sure reflection  coefficient  at  an  axial  position  x»  is  Ru^s).  The 
cross-sectional  area  of  the  duct  i between  x — Xi  and  x = x*+i 
is  Ai.  The  menu  speed  of  sound  is  c and  the  mean  density  is  p. 
At  xo,  we  assume  that  the  pressure  reflection  coefficient  J?UQ(s) 
satisfies 


|Ru0(s)l  <1  in  Real(s)  > 0 
•(«*oW)  = 0-  (75) 


x i 


Figure  11:  upstream  pipework  system 


These  assumptions  are  true  for  a choked  end  at  x = xq  for  in- 
stance. In  the  following,  we  show  by  induction  that  for  all  t, 
|RU;(s)|  < 1 in  R,eal(s)  > 0,  and  n*(Ru,(s))  = 0. 

At  x = Xi,  the  reflected  wave  fi  is  related  to  the  incoming  wave 
9i 


/,{t)  = e-2a;i4-^ftu.(s)|9i(i)  j.  (76) 

Writing  the  continuity  of  pressure  and  mass  flux  across  x = Xi+j, 
and  using  the  boundary  condition  (76)  leads  to  the  following  re- 
lationship between  the  reflected  wave  /*+ 1 and  the  incoming  wave 
</t+i  in  the  duct  i + 1: 


/i+l  (0  = e~2s  «».•+,(»)  [Si+i  (*)] . (77) 

where 


Ru,+ i(s)  — 

Ki  = 


Xi 


K,  + RUt(s)e-s^ 
1 + KiRUi(s)e~sTi 
Aj-n  — Aj 
A i + 

— 2(gj+i  - Xi) 

c 


(78) 

(79) 

(80) 


By  inspection  \Ki\  < 1.  Furthermore,  from  (78),  one  deduces 
that  n*(Rui(s))  = 0 implies  that  n*  (RU{+l(s))  = 0 after  a Pade 
expansion  [M/M]  for  e~STi  is  made.  Let  us  write  RUi  in  the  form 
RUi(s)  = |i2Ut.(s)|e"'sai  in  equation  (78).  Then  one  easily  obtains 


IP  ,-x,2  Kf  + [flM,.|2e-Mi+2iri|flu,(5)|e-^  cos(ft) 

1 Ui+lK  }]  \ + K?\RUi\2e-2*i+2Ki\RUi(s)\e-siCGs(ei)  K J 

where  9i  = hnag(s)n  + ai  and  Si  = Real(s)r{ . 

The  difference  P between  the  numerator  and  the  denominator 
of  | RUt+1  (s)|2  is  equal  to 

P = (1  - \2e~2Reatis)ri  - 1)  (82) 

Therefore,  in  Real(s)  > 0,  |RU{(s)|  < 1 implies  that  P < 0 and 
hence  |f7Ui+l  | < 0. 

Now  we  assume  that  at  x = xo,  the  high  frequency  gain 
ho  of  Ru0(s)  satisfies 


IM  < l.  (83) 

(this  is  true  for  a choked  end).  We  show  by  induction  that  for  all 
?,  the  high  frequency  gain  /i,  of  RUi(s)  satisfies  \hi\  < 1. 

FYom  (78),  after  a Pade  expansion  [M/M]  for  e~STi  is  made, 
one  deduces  that 


Ki+hi{- l)M 
1 + Kihi(—l)M 

and  hence 


(84) 


\h 


f+.i 


Kf  + h?-  2(-l)Mg;/ii 
l + 2/Ct/i.(-l)M  + Kfhf 


(85) 


The  difference  Q between  the  numerator  and  the  denominator  of 
ht+i  is 


n 
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iMT 


b j+1 1 


a J 
b j 


zj+l  zj 


: closed-loop  stability 
guaranteed 


Figure  12:  downstream  pipework  system 


Figure  13:  sign(k\ ) and  values  of  u>uTtot  for  which  a plant  of  order  2 is 
guaranteed  to  be  stabilised  by  the  fixed  controller  described  in  figure  4 

Q = (1  - *?)(/><  - !)•  (86) 

Hence,  |/i;|  < 1 implies  that  Q < 0 and  therefore  \hi+i\  < 1. 

Similarly  for  a pipework  system  downstream^  dame,  ended  Appendix  D We  start  from  equation  (70); 


by  a nozzle  (see  figure  12),  we  assume  that  at  2 = zq , the  pressure 
reflection  coefficient  Rd0  satisfies  the  assumptions  (75)  and  (83), 
which  is  true  if  the  nozzle  is  compact  [20].  By  induction,  it  can 
be  shown  that  at  any  position  z = zjt  the  downstream  pressure  where  Wm(s)  = Wcio(8)(s  4 a)  is  SPR  and  da(t)  = ^[d(f)]. 
coefficient  Rdj  and  its  high  frequency  gain  hj  satisfy 


n*(Rd.)  = 0 

| Rdj  (s)  | < 1 in  Real(s)  > 0 

\hj\  < 1* 


Pref(t)  = Wm{«)e-Jr“>!(kT(t).da(t)], 

» = Wci0(s)(s-f-a)  ii 
The  adaptive  law  is  chosen  as 

k(t)  = - Pref(t)d&(t  - not) 


(93) 


(94) 


(87) 


Appendix  B.  The  full  form  of  the  coefficients  in  equation  (7) 


X = 


+ Mf  - Mf(l  - A?i)|  (ff  - 1) 


l+j^aa. 

P2C2 

Cl  1^1  1 4 P2 


( l+Mi  (2-f)+M?(1-f) 


A?!##--1 

Ci  -y-J  1 P2 


(88) 


(89) 


(a  positive  sign  for  ko  is  assumed). 

Equation  (93)  can  be  expressed  in  a state-variable  representa- 
tion as 

x = Ax(t)  + (s  + a)  Jb(kT(f  - rtot)  da(f  - r<0t))] 

Pref(t)  - hT.x(t).  (95) 

x is  the  ‘state  vector’  of  the  system  and  A is  a matrix.  (A,  b,  h) 
is  the  ‘state  representation5  of  Wo . 

We  note  that  equation  (95)  can  be  rewritten  as 

x(f)  = Ax  A (s  4 a)  [b  (kT(t).da(t  - r*ot)}  j 


-(*  4 a) 


Appendix  C.  Consider  an  unstable  plant  of  order  2,  defined 
by  Zo(s)  = 1 and  Ro(s)  — (s  — cxu)2  4 a;2  (au  < 0 for  an  unstable 
mode),  and  whose  total  time  delay  is  rt0e-  Using  the  controller 
defined  in  figure  4,  the  denominator  of  the  closed-loop  transfer 
function  is  easily  calculated  analytically  and  is  found  to  be: 


bda  ( t - Ttot ) 


u \ 

J k(i  4-  v)dv 


(96) 


Using  equation  (94),  this  leads  to 

x(t)  = Ax4($4a)  Jb  (kT(t).da{£-not))] 


where 


Rci{s)  = $3  -f  as2  4 bs  4c 
sin  (couTtot) 


(90) 


-K«4a) 


bdaT(t-Tiot)  f J Prc/(<  + ^)da(t4^~Ttot)di> 


(97) 


a = ~2<ju  4^c  + ^2^^o^i 

OJu 

b = (j\  4wJ -2au(zc  + k2)+kokie-<T'tTtoi 
c = (al  +u}l)(zc  + k2)+kok\e~a'xTtot  ~ ^ 


As  in  Burton  [7]  and  Niculescu  [23]  [22],  the  Lyapunov  function 
candidate  is  chosen  as 


LJU 


- Sm(uJuTtot)  + COS(u>uTtot) 


sin(wii  rtot ) 4 cos(u>u  not ) 


Assuming  that  \au\  « (realistic  assumption:  the  growth 
rate  is  much  smaller  than  the  frequency  of  the  unstable  mode), 
the  Routh  criteria  applied  to  Rci(s)  shows  that  Rci{s)  is  stable  for 
some  kz  > 0 if 


k\  sm(L0uTtot)  < 0 

Zc 


(92) 


1/  tan(wuTtof)  — — < 0. 

OJu 

Consider  the  angle  0 defined  by  tan1(g)  = ^ and  sin(0)  > 0. 
Then,  as  illustrated  in  figure  13: 

• for  6 < ujuTtot  < ?r/2  and  kx  < 0,  Rci(s)  is  stable. 

• for  9 4 7T  < LJuTtot  < 27t  and  k\  > 0,  Rci(s)  is  stable. 

Hence  the  periodic  stability  bands  shown  in  figure  5a.  Also 
note  the  result  that  kx  required  for  control  satisfies  sign(ki)  — 
- sign(  sin  (u>u  rtot ) ) • 


V,  = xT(f)Px(t)  + kr(t)k(t)  + j J ||k«)||2#di/.  (98) 

(91)  “r*ot  l+u 

Using  equations  (94)  and  (97),  equation  (98)  leads  to  the  time- 
derivative 

Vt  = xT(ATP4PTA)x42xT(t)P(54a)  |fo  ^kT(t).da(t-rt0*))] 

+2xT  (t)P(s-\-a)  bda(t  — Ttot)  ^ J PTef(t  + v)da(t+V-Ttot)dL' 

0 * 

/ J(Pr.f(0)2\\dM-rtot)fdidu\  (g»; 


-2Pre/(t)kr(t).d.(e-re„e)+- 


■TtotHv 


Since  Wrn(^)  is  SPR,  lemma  2.4  [21]  can  be  used:  given  a matrix 
Q symmetric  strictly  positive,  there  exists  a matrix  P symmetric 
strictly  positive,  such  that 


22  12 


AtP+PtA  = -Q,  P(s- fa)b  = h,  (100) 

which  leads  to 


Vj  = -_xTQx_f-2Prcy(t)da(t— Tioi)  J Pref(t  Av)daL(t+v-Ttot)dvj 

o 

+ J [l|Pre/(i)da(t-Ttoe)||2~||Pre/(«+x/)da(«+ty~noi)||2]  dv.  (101) 
Denoting 

y = Pref(t)d&(t  - Ttot)>  z = Pref(t  + v)da(t  + v - rtot),  (102) 
equation  (101)  can  he  rewritten  as 

o 

Vi  = -~xtQx  + J ^2yTz  + yTy  - zTz)  dv 

— Ttot 
0 

= -xTQx-  J (||z  - y||2  - 2yTy)  Av 

— Ttot 

< —xtQx  + 2r(ot!|Pre/(t)da(t  - Tio()l|2 

< - 27Vo*||da(t  - Ttot)  ||) 

< 0 for  ||da||  < k/2rt0t-  (103) 

Therefore,  the  stability  of  the  delayed  plant  under  control  is 
guaranteed  for  all  Ttot-  A domain  of  attraction  results,  given  by 
1 | Pref  da ||  < k/2n0t]  this  domain  reaches  IRn  as  n0t  -A-  0. 
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PAPER  -22,  S.  Evesque 


Question  (M.  Mettenleiter) 

In  order  to  guarantee  the  correct  update  of  the  algorithm  you  use  a Lyapunov  function. 
How  does  this  function  depend  on  the  model  of  your  combustion  system? 

Reply 

The  Lyapunov  function  is  not  model-based  at  all.  Its  existence  relies  only  on  the  four 
assumptions  made  on  our  system  (see  paper),  but  no  other  information  on  the  system  is 
required  in  the  Lyapunov  stability  analysis. 
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